#include "KDMPrimaryGeneratorAction.hh"
#include "KDMGlobals.hh"
#include "G4Event.hh"
#include "G4PrimaryVertex.hh"
#include <G4ParticleTable.hh>
#include <G4ParticleDefinition.hh>
#include <cmath>
#include <fstream>
#include <iostream>
#include <iomanip>
#include <sstream>
using namespace std;

KDMPrimaryGeneratorAction::KDMPrimaryGeneratorAction(int photonCount) {
	rayGun = new G4ParticleGun(1);
	G4ParticleTable *particleTable = G4ParticleTable::GetParticleTable();
	G4ParticleDefinition *particle = particleTable->FindParticle("gamma");
	rayGun->SetParticleDefinition(particle);
	rayGun->SetParticleTime(0);
	numOfPhotons = photonCount;
	randEngine = new CLHEP::MTwistEngine(89347222);
	cout << endl;
}

KDMPrimaryGeneratorAction::~KDMPrimaryGeneratorAction() {
	delete rayGun;
	delete randEngine;
}

void KDMPrimaryGeneratorAction::GeneratePrimaries(G4Event* event) {
	eventNum++;

	double z = -2.75 * m + 1.6 * m * randEngine->flat();
	rayGun->SetParticlePosition(G4ThreeVector(0, 0, z));
	double polLow = atan(.141421356 * m / (-z));
	double polarAngle = (atan(.9 * m / (-z)) - polLow) * randEngine->flat() + polLow;
	double azimuthalAngle = randEngine->flat() * 6.28318531;
	rayGun->SetParticleMomentumDirection(
			G4ThreeVector(cos(azimuthalAngle) * sin(polarAngle), sin(azimuthalAngle) * sin(polarAngle), cos(polarAngle)));
	photonEnergy = (0.5 + 4.5 * randEngine->flat()) * MeV;
	rayGun->SetParticleEnergy(photonEnergy);
	rayGun->GeneratePrimaryVertex(event);

	stringstream ss;
	ss.width(15);
	ss.fill('0');
	ss << left << ((double) eventNum) / numOfPhotons * 100;
	cout << "\r" << ss.str() << " % done" << flush;
}
